Vortex arrays in a rotating superfluid Fermi gas 
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The behavior of a dilute two-component superfluid Fermi gas subjected to rotation is investigated 
within the context of a weak-coupling BCS theory. The microscopic properties at finite temperature 
are obtained by iterating the Bogoliubov-de Gennes equations to self-consistency. In the model, alkali 
atoms are strongly confined in quasi-two-dimensional traps produced by a deep one-dimensional 
optical lattice. The lattice depth significantly enhances the critical transition temperature and the 
critical rotation frequency at which the superfluidity ceases. As the rotation frequency increases, 
the triangular vortex arrays become increasingly irregular, indicating a quantum melting transition. 
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Great experimental strides have been made over the 
past year in the goal to form a BCS-like superfluid 
state with ultracold neutral Fermi gases. Through the 
use of magnetic field-induced Feshbach resonances jjj, 
0, El, strong interactions between quantum degenerate 
fermions in two different hyperfine states have been in- 
duced 0, IE IE 0- Molecules result when the result- 
ing interactions are repulsive || 0, El El E3 ! though 
Cooper pairing at high temperatures (within an order of 
magnitude of the degeneracy temperature Tp) is widely 
expected for attractive interactions near the Feshbach 
resonance lr3. IT3. t he state actually produced in recent 
experiments 

en 

16| remains to be fully elucidated. 

Another proposed approach to high-temperature su- 
perfluidity in these systems is to confine the gases in op- 
tical lattices El El • The flattening of the energy bands 
as the lattice depth increases ensures that a large num- 
ber of atoms participate in the pairing even for relatively 
weak interactions. A deep one-dimensional (ID) lattice 
corresponds to an array of quasi-2D traps, whose large 
effective axial confinement can enhance the interaction 
strength sufficiently El to raise the superfluid transition 
temperature T c to within 10% of T F 20]. 

Rotating the superfluid around the optical lattice axis 
should yield vortex arrays akin to those observed re- 
cently in trapped Bose- Einstein condensates |2l). Few 
clear signatures of superfluidity exist for dilute Fermi 
gases J22l: th e first evidence has been obtained only very 
recently [23I E3 . l2fij . The experimental observation of 
quantized vortices would be a 'smoking gun' for the pres- 
ence of superfluidity in the system, though significant de- 
pletion of particles in the vortex cores is expected only 
in the strong-coupling limit j2(|. Vortices in the array 
of weakly coupled two-dimensional (2D) traps serve as 
a clean model of the pancake vortices in strongly lay- 
ered high-Tf, superconductors such as r^S^CaC^Os+a 
(BSCCO) [23. Furthermore, because a small fraction of 
the total number of atoms reside in any individual lat- 
tice site, one might expect to readily access the quantum 
Hall limit at high rotation frequencies where the num- 
ber of vortices becomes comparable to the number of 



atoms HHH. 

To my knowledge, the calculations presented below 
constitute the first microscopic determination of vortex 
arrays in inhomogeneous Fermi superfluids. The central 
results of the work are: (1) the critical rotation frequency 
for the cessation of superfluidity (the rotational analogue 
of the critical velocity v c in superconductors) scales with 
T c ; (2) the vortex cores for superfluids in the BCS limit 
that are confined in ID optical lattices should be directly 
visible experimentally at low temperatures; and (3) the 
vortex lattice is expected to melt at zero temperature at 
high rotation frequencies approaching the expected quan- 
tum Hall transition. 

The calculations were inspired by experiments with at- 
tractive Fermi gases confined in optical traps 0, El ■ The 
experiments were performed at magnetic fields above but 
near the Feshbach resonance at 860 G |30j at which the 
s-wave scattering length between atoms in the hyper- 
fine states IF = |, mjr = ± |) (labeled a =T,J, below) 

is a 10^a (a = 0.0529 nm). With iV = 7.5 x 10 4 

atoms in each species Q, the Thomas- Fermi (TF) ap- 
proximation N = ^^f-(^f-) 3 can be inverted to yield the 
TF chemical potential /i^F ~ 25hu p and the Fermi mo- 
mentum at the trap center hk F sw ^/2to^tf; one then 
obtains k F \a\ 3> 1 which implies very strong coupling. 

In contrast, the BCS mean-field approximation is ap- 
plicable only for k F \a\ < 0.6, beyond which one needs 
to include pairing fluctuations 31] . The BCS transi- 
tion temperature is given approximately by the uniform 
3D expression [12 T c = {8e^~ 2 /tt)h tf exp[-7r/(2fc£|a|)], 
where 7 = 0.5771 • • • is Euler's constant; this does not 
take into account polarization corrections (not considered 
in the present work), which lower T c by a factor of ap- 
proximately 2.2 j33j - Decreasing |o| such that k F \a\ = 0.6 
yields T c ss lAhujp ~ 0.06Tf which remains too low to 
easily access experimentally. 

This value for T c can be substantially increased if the 
laser were retro-reflected to form a ID optical lattice, 
with potential Vf a t = tiEr sin 2 (27rz/A) where n is the 
lattice depth in units of the recoil energy Er = ^(x) 2 - 
Near the center the potential may be approximated as 
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|s) where d z = y/h/mLj z , not in- 



quadratic in z, with the ratio of the effective axial and 
radial frequencies a = uj z /u p = ^/n(2ird p / A) 2 = 20 — 80 
where d p = y^h/muip and the lower and up per lim- 
its correspond to A = 1064 nm, n sa OAEr [iff and 
A = 10.6 ^m, n w 50 000£^ 0. Neglecting inter-well 
tunneling, all of the atoms occupy the axial oscillator 
ground state in the quasi-2D limit fi <C hu) z . A straight- 
forward calculation [34| gives the superfluid transition 
temperature for a uniform quasi-2D system ksT 2T> = 
(2eT/7r)/x2Dexp(- 
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eluding polarization effects which reduce T C 2D by a fac- 
tor of e |2jj. In the calculations below, /1 = lOhujp, 
dp = 1 /im, a = — 2160ao (corresponding to the triplet 
scattering length that would be obtained for an external 
bias field around 2000 G far above the Feshbach reso- 
nance), and 20 < a < 40. Choosing a — 40, one obtains 
T 2D 0.16 /iK« T F /5. This significantly improves T c 
while fcf/^laj « 0.5 remains in the weak-coupling limit. 

The calculations are based on a BCS mean-field ap- 
proximation to the full 2D interaction Hamiltonian 35], 
where the particle density and gap functions are de- 
fined by the thermal averages ^(x) = (x)^ CT (x)) and 
A(x) = — ^'(V'f (x)^j(x)), respectively (here and below 
x = (x,y) corresponds to a 2D vector). The ultraviolet 
divergence in the definition of the superfluid gap is reg- 
ularized using the pseudopotential method |36L |3J| , giv- 
ing rise to a regularized coupling constant g' . An equal 
population of the two hyperfine components — is 
chosen to maximize T c |32j ] . Diagonalizing the mean-field 
Hamiltonian in a frame rotating around the lattice axis 
at angular frequency fi, one obtains the Bogoliubov-de 
Gennes (BdG) equations [3^] 
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n-n A(x) 

A*(x) -(H-fi)* 



u„(x) 



E„ 



u„(x) 



(1) 



Here H = ~^V 2 + \mu? p (x 2 + y 2 ) + g'n a {x) - QL Z , 
where L z — ih(yd x — xd y ). The density and gap func- 
tions are n a (x) = [\u n \ 2 f(E n ) + \v n \ 2 (l - /(£„))] 
and A(x) = -g' J2 n u nV*(^ " V{ E n)), respectively, and 
must be iterated to self-consistency together with the 
quasiparticle amplitudes u n and v n appearing in Eqs. 
to obtain equilibrium solutions; the sums run over posi- 
tive E n and f(E n ) = (e E "/ kBT + l)" 1 . 

The BdG matrix was evaluated numerically using a 
DVR based on Hermite polynomials with up to 100 func- 
tions/points each in x and y. The eigenvectors for en- 
ergies up to E n (max) were obtained using the routine 
pzheevx on a 32-processor Xeon cluster. For the largest 
grid, each diagonalization took approximately 50 min- 
utes. The procedure was deemed converged when the 
magnitudes of A and n CT at successive iterations were 
smaller than some predefined tolerance; for large fi as 
many as 3000 iterations were required to verify that the 
ground state had been obtained. The sums over E n were 
found to fully converge for _E(max) ~ 3[i. N a was on the 
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FIG. 1: (T, fi) phase diagram for the quasi-2D system as a 
function of a = Cj z /lo p . Regions to the left (right) of a given 
curve correspond to superfluid (normal) states. Note that 
Tiujp/kB w 80 nK, so that T c (fi = 0) fa 200 nK for a = 40. 



order of 3000 for fi = 0.98w p for a = 40. 

Because the effective coupling strength g' oc \fa oc n 1 / 4 
(where n is the lattice depth in recoils) is a parameter 
adjustable by varying the trapping laser intensity, it is 
useful to explore the (T, fi) superfluid phase diagram for 
various a, subject to the constraints fi/huip < a (single 
axial mode approximation) and k F \a\ < 1 (weak cou- 
pling). Results for a = {20, 30, 40} are shown in Fig. ^ 
the superfluid transition was obtained assuming cylindri- 
cal symmetry (when A = there are no vortices) with 
random points checked using the full 2D code. The nu- 
merical value of T c for fi = is comparable to the TF 
value of T 2D given above, in spite of the small total num- 
ber of atoms and the inhomogeneous density. Superflu- 
idity also ceases at zero temperature at a critical angu- 
lar frequency fi c . This is entirely due to the breaking 
of Cooper pairs by the quasiparticles' rotational motion, 
and is the neutral analog of the critical velocity v c in uni- 
form superconductors subjected to a voltage drop |35|. In 
the latter case, v c exec A(T = 0) oc T c ; the numerical re- 
sults clearly indicate that while the critical frequency is 
proportional to T c , it is much smaller than u p unless T c 
is enhanced by the trap anisotropy. 

The immediate question is: can fi c be made sufficiently 
large that one or more vortices can be stabilized? If fi is 
too small, the vortex-free state is energetically favored, 
and for fi ~ fi c the magnitude and spatial extent of A 
are too small to support vortices whose characteristic size 
is the local healing length £ = H 2 kF /rrnrA. No vortices 
were observed for a = 20, while a = 30 yielded a small 
number, fewer than 10. As shown in Fig. [21 the a = 40 



3 



□□□□ 

nn 

:::: 

FIG. 2: Density plots for the in-trap gap function A(x) (a- 
c) and number density n,x(x) (d-f) for a = 40 and T = 0. 
Figures (a,d), (b,e), and (c,f) correspond to Q./u> p = 0.4, 0.5, 
and 0.7, respectively. Each box is 16 /xm on a side. 




FIG. 3: Density plots for the in-trap gap function A(x) (a,b) 
and number density n a (x.) (c,d) for a = 40 and Q./u) p = 0.4, 
for which kBT c /huj p ~ 2.03. Figures (a,c) and (b,d) corre- 
spond to ItBT/hujp = 1 and 2, respectively [the T = cases 
are shown in Fig|!Ja,d)]. Each box is 16 /im on a side. 



case, with ft c 



is able to support large numbers. 



For small ft, the vortices distribute themselves into a 
regular triangular pattern; however, for increasing an- 
gular frequency the arrays become disordered, so that 
by ft = 0.7ujp the T = lattice has completely melted. 
These results were independent of the initial guesses for 
A and n a and of the convergence criterion. The vortices 
tend to be found on circles defined by the boundary be- 
tween orbitals of angular momentum m and m + 1, and 
concomitantly exhibit considerable azimuthal distortion 
(the effect is most noticeable in the number density). 

The zero-temperature melting of a 2D vortex array 
was suggested some time ago in the context of the lay- 
ered high-T c superconductors without impurities |3SL l39| . 
According to the Lindemann criterion, melting occurs 
when the zero-point amplitude of vortex position fluc- 
tuations is some critical fraction cl of the intervortex 
separation £ v ; for most materials, ci « 0.1 - 0.2. The 
quasiparticle bound states in the vortex core account for 



the vortex fluctuations, and have a spatial extent ~ £. 
The vortex separation £ v for a given ft is determined by 
matching the angular velocity of the normal-state atoms 
with the average circulation k — h/2m of the superfiuid 
vortices, which yields the vortex density n v — 2ft/fc or 
£ v = rs 2d p = 2 /im which depends only weakly 

on ft. For ft = {0.4, 0.5, 0.7} shown in Fig. H^a-c), the 
numerics yield £/£ v {0.15, 0.17, 0.2}, consistent with 
the Lindemann criterion. Indeed, visual inspection of 
the lowest-energy quasiparticle amplitudes u n and v n 
reveals that they are highly localized near the vortex 
cores for ft = 0.4 but begin to overlap those of adja- 
cent vortices for larger ft. Furthermore, for systems in 
the lowest Landau level the relation N/N v ~ c7/ is ex- 
pected (2^; f° r the values of ft considered above one 
obtains N a /N v = {12.6, 8.7, 5.5}. An important ques- 
tion beyond the scope of the present work is to elucidate 
the relationship (if any) between the vortex array melt- 
ing observed here and any impending transition into a 
quantum Hall state (not accessible with the current for- 
malism). Fractions N/N v < 6 have been predicted to 
favor a quantum Hall state in Bose gases |28| . 

Perhaps the most important result of this work is that 
the vortex cores are faint but nevertheless clearly vis- 
ible as depressions in the particle density, particularly 
for lower ft and T. For ft = 0.4 and 0.7, the density 
in the vortex core is reduced by approximately 1/2 and 
1 /4, respectively. Such large depressions were previously 



4 



thought only to occur far from the weak-coupling BCS 
limit relevant to these calculations |3^ . 

As shown in Fig. [3] the results for a given f2 are 
strongly affected by the temperature. As T increases, 
only atoms in the vicinity of the trap center participate 
in the pairing because the local value of T c is lowest at 
the low-density surface of the cloud; alternatively, most 
of the pair-breaking quasiparticle excitations occur at 
the surface where the effective potential is lowest and 
the tangential particle velocities are largest. Similar re- 
sults have been obtained with a local-density Ginzburg- 
Landau theory (ijj. Because the healing length, which 
governs the vortex core size, also diverges as the temper- 
ature increases, fewer vortices can fit into the reduced 
superfluid region. The visibility of vortices in the par- 
ticle density is reduced at finite temperature due to the 
thermal occupation of core states. 

The most important outstanding issue to be addressed 
in future work is the vortex coherence between the wells 
of the optical lattice. For deep lattices such as are consid- 
ered here the Josephson coupling between sites becomes 
small |4l|. In this regime, vortex lines could break up 
into disconnected 'pancake' vortices in each layer, ren- 
dering them unobservable experimentally. In the high-T c 
superconductors, this vortex decoupling occurs when the 
anisotropy parameter 7 = \Jm z /m p > 150 |2tI [4^ (m z 
and m p are the effective masses perpendicular and par- 
allel to the layers), though a full description of such a 
decoupling transition is not yet complete . Assuming 
m z = h 2 1 (9 2 £fc/9fc 2 ) where £k is the dispersion of the 
lowest-lying band near the trap center, the decoupling 
7 corresponds to an optical lattice depth approaching 
7QEr. This value is much lower and higher than those of 
Refs. M and flil respectively, suggesting that the melting 
of well-defined vortex lines should be directly observable 
in experiments employing Nd:YAG (but not CO2) lasers. 
A full calculation incorporating interwell tunneling re- 
mains necessary; the disappearance of vortices observed 
experimentally could be due not only to vortex decou- 
pling but also to the onset of a quantum Hall phase, 
which is expected to restore the cylindrical symmetry 
broken by the vortices. 
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